* Program to compute ASL of difference of means using naive moving block bootstrap
* and HAC-corrected standard errors

* This program: 1980-1992 vs. 1994-2006 samples, RMSE-based NSRs


calendar 1947 1 4
allocate 2020:04+5000


open data table1_reg_test_data.rat
data(format=rats) / gb_revgapsq_f13



seed 43892



* Set bootstrap and related parameters

compute ndraws=5000
compute nblocks=4
compute haclags=4
set t_orph_rmse = %NA


set dum94q1 = 0.0
set dum94q1 1994:01 2006:04 = 1.0

linreg(robusterrors,lags=haclags,lwindow=neweywest) gb_revgapsq_f13
# constant dum94q1


statistics(noprint) gb_revgapsq_f13 1980:01 1992:04
compute mean1=%mean
statistics(noprint) gb_revgapsq_f13 1994:01 2006:04
compute mean2=%mean
statistics(noprint) gb_revgapsq_f13 1980:01 2006:04
compute meanfull=%mean

disp "Difference in means = " mean2-mean1


set zzz 1980:01 1992:04 = gb_revgapsq_f13-mean1+meanfull
set yyy 1994:01 2006:04 = gb_revgapsq_f13-mean2+meanfull


do iter=1,ndraws

  set b_zzz = %NA
  set b_yyy = %NA

  boot(replace,block=nblocks,method=overlap) b_entry1 1980:01 1992:04
  set b_zzz = zzz(b_entry1)

  boot(replace,block=nblocks,method=overlap) b_entry2 1994:01 2006:04
  set b_yyy = yyy(b_entry2)

  linreg(robusterrors,lags=haclags,lwindow=neweywest,noprint) b_zzz
  # constant
  compute mean_bzzz=%beta(1)
  compute se_bzzz=%stderrs(1)

  linreg(robusterrors,lags=haclags,lwindow=neweywest,noprint) b_yyy
  # constant
  compute mean_byyy=%beta(1)
  compute se_byyy=%stderrs(1)

  set t_orph_rmse iter iter = (mean_byyy-mean_bzzz)/sqrt((se_bzzz**2.0)+(se_byyy**2.0))


end do iter



* Tabulate bootstrapped t-statistics

statistics(fract) t_orph_rmse


